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Microcanonical critical droplet theory and molecular dynamics simulations are used to examine 
static coexistence between solid and liquid phases in nanoscale lead clusters. It is shown that the 
theory predicts the existence of a metastable coexisting state above a critical cluster radius R± , with 
this state becoming stable for clusters of radius R > R2. Molecular dynamics simulations of lead 
clusters confirm the existence of stable coexisting states in 1427-atom and 2057-atom clusters but 
find no stable coexisting state in a 931-atom cluster. 
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I. INTRODUCTION 

The coexistence of solid and liquid phases in atomic 
clusters has been of experimentally and theoretical^^ 
interest for some time. In large clusters, as in bulk mate- 
rials, static coexistence is expected theoretically and is 
readily observed in both experiments^ and simulation^. 
However in sufficiently small clusters, simulations find 
that coexistence is dynamic i.e. clusters fluctuate in time 
between fully solid and fully liquid states^. Small clus- 
ters avoid static coexistence because of the prohibitive 
cost of forming an interface. This causes an S-bend in 
the microcanonical caloric curved and the corresponding 
negative heat capacities which been observed in small 
sodium clusters^.. 

However, it is currently not known at what size static 
coexistence gives way to dynamic coexistence. The time- 
scale on which dynamic coexistence occurs will be gov- 
erned by the height of the energy barrier separating the 
solid and liquid states of the cluster—. However, if a coex- 
isting state exists, then this coexisting state should also 
be dynamically accessible. This is not seen in at least 
some small cluster o 3 : 10 . Evidently, in sufficiently small 
clusters coexistence is unstable. At what sizes do states 
of static coexistence occur? 

To determine when phase coexistence will be stable, 
one needs to know the relative cost of forming an in- 
terface between phases in a given cluster. A common 
approach to studying coexistence in clusters is to apply 
critical droplet theory in the capillarity approximation 
(i.e. where cluster material parameters, especially sur- 
face energies, are approximated by their bulk values). 
For instance Reiss et all used the capillarity approxi- 
mation to develop a critical droplet theory of melting in 
the canonical ensemble, with the (spherically-symmetric) 
solid fraction of a cluster acting as an order parameter. 
Nielsen et al& developed a similar critical droplet theory 
in the microcanonical ensemble. They demonstrated that 
a stable static coexisting state can exist in this theory 
near the melting point, and validated this with constant- 
energy simulations of a large copper cluster. Cleveland et 
ali developed an aspherical model which dealt with static 
coexistence where the liquid only partially wets the solid. 



Thus, in principle, critical droplet theory can be used to 
discuss coexistence in a variety of circumstances. 

To study coexistence we will work in the microcanon- 
ical ensemble as constant temperature ensembles tend 
to suppress the inherent inhomogencity in coexisting 
systems-' 1 '. Constant energy conditions are realised in 
isolated systems with poor thermal coupling to their en- 
vironment, such as clusters in an inert gas atmosphere. 
Consequently we will use microcanonical critical droplet 
(MCD) theory to examine the stability of static solid- 
liquid coexistence in small metal clusters. We will then 
use constant-energy molecular dynamics simulations to 
test the predictions of the theory. 

II. MICROCANONICAL CRITICAL DROPLET 
THEORY 

Here we use microcanonical critical droplet (MCD) 
theorj£ to examine the stability of static coexistence. 
The spherical coexisting solid-liquid cluster of radius R, 
is considered to consist of a solid core, of radius R s , and 
a liquid outer shell. The entropy (per atom), S m , of a co- 
existing solid- liquid cluster of solid fraction r\ = (R s /R) 3 
at energy (per atom) e is approximated by 

S m (e)=r,S s (e s ) + (l-r,)S l (e l ) (1) 

where S s (e s ) = S° + clog + l) and Si(e t ) = 

Sf + clog ( e ' e ^ e ' + . Here T c is the bulk melting tem- 
perature, c is the heat capacity and quantities with a 
superscript c indicate values at the bulk melting point. 
The solid energy, e s , and liquid energy, e;, are related to 
the total energy, e, for a cluster of radius R as follows: 

e = T]e s + (1 - rj)ei (2) 

+ ^ (7/, + (1st + A7ex P (-2(1 - ^R/d))^) 

where A7 = j sv — j s i — ji v (-y sv is the solid- vapour in- 
tcrfacial energy density, 7 S ; is the solid-liquid interfacial 
energy density and ji v is the liquid- vapour interfacial en- 
ergy density), p is the bulk material density and £ is 
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a correlation length which characterizes the short-range 
interaction between the solid-liquid and liquid-vapour 
interface^. The inclusion of this final term in equation 
(J2J ensures that e — + e s + 3j sv /pR as r\ — ► 1. 

We now show that MCD theory predicts that there is 
a critical cluster radius below which coexistence becomes 
unstable. For a given total energy per atom e, we can 



regard e; as a function of e s and 77: e; = ei(e s ,r]). We 
may then extremise the entropy (given by of the 
mixed cluster with respect to e s and 77 at fixed e. It is 
then straightforward show that the entropy is extremiscd 
when the following conditions hold: 1) 7} = T Sl which is 
the normal condition for thermal equilibrium, and 2) 77 
must satisfy the following equation^: 



114/3 -^L ( 7si + ^- 2(lW ' 3)m ) +( e - e ^7| +! ^ e " (1 ^ 1/3)it/5 ) ! T ! + ^ 0< + A7 e - 2(1 - 1/3) ^) = 

9 9 9 9 _ m. 



where L = e\ — e c s is the latent heat of fusion per atom. 
For sufficiently negative values of e, there are two positive 
solutions to ?/i (e) > 772(e) > 0. There is always a 
range of energies where < 772(e) < 1, as when e — ► —00, 
772(e) — > 0. In fact 772(e) represents a local minima in the 
entropy, corresponding to a barrier separating the solid 
and liquid phases. Similarly, for values of e where 771(e) 
exists and is less than one, 771(e) is a local maxima in the 
entropy, corresponding to a coexisting state (note that 
ik<0andf| J)2 >0). 

A range of energies where < 771(e) < 1 exists only 
when R > Ri, where R\ is given by: 

Ri = ^(-A+V# + B), (4) 

where A = g(^-2^(l-^)) and B = 

2(1 + {j sv - j lv ). That is if R > Ri, there is a 

range of energies where a locally stable coexisting state 
exists. At such an energy e, the coexisting state, with 
solid fraction 771(e), is separated from the liquid state by 
an energy barrier at 772(e). 

However, the coexisting state with solid fraction 771(e) 
is not necessarily globally stable. If the coexisting state 
exists at an energy e, then S m {e) > S s (e). However, 
there exist a range of energies where S m (e) > Si(e) only 
if R > R2, where 

R = 7s, - liv I 3£ - 2cr c (1 - e - L / cT <) \ 
PL [ L - cTc (l-e-^){l-^))- 

(5) 

Thus for clusters with Ri < R < R 2 , at energies where 
the coexisting state exists, the state is always metastable 
(S m (e) < 5/(e)). For clusters with R < R\, no coexisting 
state exists at any energy. 

We note that some of these parameters may be ill- 
defined for small clusters. Nonetheless the values of R\ 
and i? 2 may be calculated using the capillarity approx- 
imation i.e. using bulk values (if available) for the ma- 
terial parameters. For leadi^ we find that Ri — 0.2 nm 



(or approximately 1 atom) and R2 = 0.6 nm (or approxi- 
mately 100 atoms). This suggests that coexistence ought 
to be stable down to quite small cluster sizes, although it 
is well-known that the melting behaviour of sub- 100 atom 
clusters is quite erratic due to electron-shell effects 13 that 
are certainly not captured in the MCD model. 



III. SIMULATED CALORIC CURVES 

Using molecular dynamics simulations of Pb clusters 
we have found that coexistence becomes unstable at ap- 
proximately 1000-atoms. We simulated the caloric curves 
for three clusters: a 2073-atom icosahedron, a 1427-atom 
icosahedron and a 931-atom icosahedron. Here we use 
surface-reconstructed icosahedra which are thought to 
be stable in Pb cluster a 14 i 15 at these sizes using a glue 
potential^. At each energy the cluster was equilibrated 
for 0.6 ns, then the kinetic energy was averaged over 
a further 0.6 ns. The energy increment used was 0.2 
meV/atom with energies adjusted between constant en- 
ergy simulations by a uniform scaling of the kinetic en- 
ergy. Figure n shows the caloric curves for the 2073 and 
1427-atom clusters moving from the solid state (on left) 
to the liquid state (on right). Both clusters exhibit a co- 
existing solid-liquid state separating the solid and liquid 
phases. The curves in figure^resemble those constructed 
in Ref— with their layer-by-layer model of cluster melt- 
ing. 

Note the change in temperature at the transition from 
fully solid to coexistence and in the transition from co- 
existence to fully liquid. In MCD theory, while the tran- 
sition from coexistence to liquid is first-order, the tran- 
sition from solid to coexistence is continuous. To check 
whether the 1427 and 2057-atom clusters are superheated 
when they undergo the transition from solid to the co- 
existing state, we cooled a 1427-atom cluster from the 
coexistence region down until the cluster froze. Figure |3 
shows the reverse caloric curve (constructed using equi- 
libration times of 4 ns with an energy increment of -0.2 
meV/atom). The freezing transition takes place to well 
within 5 meV/atom of the transition from solid to coexis- 
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FIG. 1: Caloric curves for (a) a 1427-atom Pb cluster and 
(b) a 2073-atom cluster. The 1427 and 2073-atom clusters 
exhibit coexistence for a range of energies. 



tence, suggesting that the simulated caloric curve is close 
to the equilibrium caloric curve and that the transition 
at -1.715 eV/atom is a first-order transition which ex- 
hibits little rounding due to finite-size effects. The origin 
of this discontinuity at this transition may correspond 
to the creation of a "minimum" volume of liquid and as 
such the discontinuity is likely to disappear as R — > oo. 

To characterise the coexisting state, we follow Cleve- 
land et afi in using the bimodality of the distribution of 
diffusion coefficients to distinguish solid and liquid atoms. 
A sequence of snapshots at different total energies of the 
f427-atom cluster is shown in figure |3| The total energy 



FIG. 2: Caloric curves for a 1427-atom cluster produced first 
by heating (solid line), then by cooling (dashed line). 



of the atoms in each snapshot decreases from top left to 
top right to bottom left to bottom right. Each atom in 
the snapshot has been coloured either white or grey de- 
pending on the whether the diffusion coefficient of that 
atom (calculated over a carefully chosen timescale) is low 
mobility (solid) or high mobility (liquid). The snapshots 
reveal that the molten lead does not fully wet the solid. 
They also reveal that the liquid fraction of the coexisting 
cluster increases as the total energy increases. Note that 
the melt does not wet the solid as assumed in equation [21 
whereas in bulk lead, the melt wets the solid (since in 
bulk lead A7 > 0). This is a clear indication that the 
capillarity approximation has broken down at these clus- 
ter sizes. 

Figure 01 shows the caloric curve for a 931-atom clus- 
ter constructed in the same way as those for the larger 
clusters in figure ^ Note that the curve does not show a 
stable solid-liquid coexisting state but that large fluctua- 
tions in temperature occur near the melting transition in 
the 931-atom caloric curve. Using the same method for 
identifying solid and liquid regions as was used in figure|3 
one can verify that the fluctuations in temperature are 
correlated with the appearance of a liquid region in the 
otherwise solid cluster. This is due to the dynamic co- 
existence between the solid state and a coexisting state. 

Which state is globally stable? We have resolved this 
with a very long duration run near the melting transi- 
tion. Figure shows a history of the temperature (as de- 
termined using the method above) for a 931-atom cluster 
at a total energy of E = —1.688 eV/atom over approx- 
imately 15 ns. The final state of the cluster is liquid 
(achieved after approximately 8 ns), which suggests that 
the fluctuations (between and 8 ns in figure and be- 
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FIG. 3: A sequence of snapshots showing the coexisting 
solid-liquid state as the energy is increased from E = —1.712 
eV/atom (top left) to E — —1.706 eV/atom (bottom right). 
Atoms are classified by their mobilities: the dark grey atoms 
are liquid (high mobility) while the light grey atoms are solid 
(low mobility). 

tween —1.690 and —1.686 eV/atom in figure^} are non- 
equilibrium precursors to melting rather than true dy- 
namic coexistence. MCD theory certainly suggests that 
at sizes R < R2 the coexisting state is metastable with 
respect to the liquid rather than the solid. However, 
the caloric curve in figure 0] suggests that the coexisting 
state is metastable with respect to both the solid and 
liquid states, as the cluster appears to spend more time 
in the solid state. Nonetheless, we conclude that there is 
no stable coexisting state for the 931-atom cluster at any 
energy. In Refi&, we constructed a less detailed caloric 
curve for 1130-atom lead clusters by coalescence. This 
curve also showed no evidence for stable coexistence. 



IV. DISCUSSION 

The molecular dynamics simulations and MCD theory 
in the capillary approximation do not agree on the size 
at which coexistence becomes metastable in Pb clusters. 
However, the snapshots in figure [3] strongly suggest that 
the capillarity approximation has broken down at these 
cluster sizes (at least when it comes to reproducing the 
behaviour of the simulated clusters). This should not be 
surprising as atomic clusters are interesting precisely be- 
cause their properties depend strongly on size. However, 
if we are to take MCD theory beyond the capillarity ap- 
proximation, we would require a detailed understanding 
of how surface energies, and other quantities such as la- 
tent heat, depend on surface curvature and/or size. It is 
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FIG. 4: Caloric curves for a 931-atom Pb cluster. The 931- 
atom cluster does not exhibit a stable coexisting phase. The 
large fluctuations near the melting point signal the appear- 
ance of precritical liquid nuclei. The true melting point energy 
probably lies closer to the dotted line. 
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FIG. 5: Time evolution of the kinetic temperature in a 931- 
atom cluster at an energy E — —1.688 eV/atom. The white 
line shows a time averaged trace of the temperature. The 
cluster starts as solid (0-5 ns), becomes coexisting (5-8 ns) 
and finally melts fully (after 8 ns). 

not simply a matter of inserting more accurate material 
properties into JSJ. 

It is also important to recognise that the spherically 
symmetric geometry assumed in equation is inap- 
propriate for situations involving partial wetting where 
A7 < 0. In Refi, a critical droplet model was developed 
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for the partial wetting case, although it was applied in the 
canonical ensemble. Evidently a generalisation of MCD 
theory to the partially wetting case would also be desir- 
able. Similarly it would be interesting to determine by 
simulation and experiment the size at which A7 = in 
Pb clusters i.e. the size below which complete wetting 
gives way to partial wetting. Partial and complete wet- 
ting has been observed recently in isolated metal alloy 
clusters using transmission electron microscopy^. 

Finally, in the simulations the transition from solid to 
solid-liquid coexistence was discontinuous. In MCD the- 
ory this transition is continous as the energy of the co- 
existing cluster, as given by equation (J2J, approaches the 
energy of an entirely solid cluster as r\ — > 1. This limit is 
an artificial feature of equation (J2J); if this constraint were 
to be relaxed so that the energy of the coexisting clus- 
ter approached some other value as 7/ — > 1 (for instance, 
e s + 3(7;^ + j s i)/pR) then the transition from solid to 
solid-liquid would be discontinuous. Provided the energy 
of the solid cluster and that of the limiting solid-liquid 
cluster differed only by surface energy, this discontinuity 
would disappear as R — > 00. 

V. CONCLUSION 

In summary, MCD theory predicts that the static co- 
existence of solid and liquid phases can only occur in 



clusters above a certain size given by equation J5j|. Be- 
low this critical size coexistence is metastable or unsta- 
ble. Molecular dynamics simulations of Pb clusters using 
a glue potential appear to corroborate this scenario, al- 
though the transition occurs at a larger size, N2, (where 
1100 < N 2 < 1400) than that predicted by MCD theory 
using bulk material parameters in the capillarity approxi- 
mation. Indeed this discrepancy between theory and sim- 
ulation is probably due to the breakdown of this approx- 
imation at small cluster sizes. We note that Pb has an 
unusually large value of the prefactor x — (j sv — ji v )/pL 
that appears in equation JSJ due to a relatively small la- 
tent heat of fusion. Other metals with small latent heats 
(which are correlated with low melting points), such as 
Hg or the alkali metals, will also have large values of x. 
It may be that in such metals coexistence is forbidden 
at quite large sizes. In most metals, however, with rela- 
tively large latent heats of fusion, we might expect static 
coexistence to remain stable down to rather small sizes. 
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